High species diversity of Phintella and Phintella‐like spiders (Araneae: Salticidae) in Vietnam revealed by DNA‐based species delimitation analyses

Abstract Salticidae (jumping spiders) usually exhibit pronounced sexual dimorphism in adult morphology, particularly body coloration and size and shape of the first legs. Consequently, the male and female from the same species might be erroneously assigned to different species or even different genera, which could generate synonymies in classification if only morphological data were used. Phintella is a species‐rich genus of Salticidae, which currently exhibits 76 named species. However, the male–female counterpart is unknown for nearly half of the species. In this study, we used a molecular approach to delineate the species boundaries for Phintella and Phintella‐like specimens collected in Vietnam, using morphological information as supporting data. We used three gene fragments (mitochondrial COI, 16S‐ND1, and nuclear 28S) and biogeographical considerations for species delimitation. A total of 22 putative species were recognized: 18 species of the genus Phintella, one species of the genus Lechia (L. squamata), and three species of the genus Phinteloides. Eleven undescribed species were discovered, of which seven have a male–female combination, two species have only males, and two species have only females. The crown age of Phintella was estimated at the Serravallian stage of the Miocene after the increase of species number around 16 MYA. The crown ages of most putative species recognized in this study were estimated in the Pleistocene, and the divergence among sister species likely occurred from the mid‐Miocene to the Pliocene. Our ancestral range reconstruction results showed that the diversification of our ingroup was governed by progressive dispersal events, i.e., Phintella and their related species in Vietnam diversified while expanding their range on the continent. Our results provide fundamental biodiversity data for a high‐diversity genus in Vietnamese Phintella spiders.


| INTRODUC TI ON
The order Araneae (spiders) is a major group of predatory arthropods in terrestrial ecosystems (Foelix, 1996).The global spider community consumes 400-800 million tons of prey annually, approximately 0.1% of the global terrestrial net primary production (Nyffeler & Birkhofer, 2017).Salticidae, the most species-rich spider family (Bodner & Maddison, 2012;Richman, 1982), consists of 6640 described species and 681 genera, holding 12.8% of all species of Araneae (51,917 species, World Spider Catalog, updated February 15, 2024).In spiders, the morphology of male and female copulatory organs (male palp and female epigyne) has features related to mechanical reproductive isolation, which have been intensively used in delineating the boundaries among species or higher taxa (Foelix, 1996).However, the current classification of Salticidae, especially the records in Southeast Asia, is primarily based on morphology, which may involve many errors and obscurities in taxon recognition (e.g., Phung et al., 2016).As Salticidae usually exhibit distinct sexual dimorphism in adult morphology, particularly body coloration, body size, and shape of the first legs, conspecific males and females might often be erroneously assigned to different species or genera.Given that genetics or behaviors are not available for matching conspecific males and females, many synonymies might be hidden in the current classification (Logunov & Jaeger, 2015;Phung et al., 2016;Żabka, 1985).
Species delimitation using a combination of molecular phylogenetic analyses and morphological examination has been widely applied to various animal taxa (Derkarabetian & Hedin, 2014;Diaz et al., 2015;Gonzalez et al., 2014;Kergoat et al., 2015;Montagna et al., 2016;Su et al., 2016).Arachnologists are vigorously elucidating the species diversity of spiders and updating the conventional classification using molecular data while taking morphological identification into consideration (Ballarin & Eguchi, 2023;Li et al., 2022;Lo et al., 2021;Macharoenboon et al., 2021).Such approaches are practical in solving the taxonomic problems and obscurities in the species recognition and classification of Salticidae, which are primarily due to remarkable sexual dimorphism in adult morphologies (Bopearachchi et al., 2022;Kanesharatnam & Benjamin, 2021;Phung et al., 2016).Moreover, Maddison et al. (2020) demonstrated the effectiveness of a genomic approach to draw species boundaries in Salticidae with an example of the tribe Baviini, whereas the genus-level definition and species-level recognition are still primarily based on morphology (Prószyński, 2016a(Prószyński, , 2016b(Prószyński, , 2019)).
The jumping spiders in the genus Phintella are no exception to the taxonomic obscurities originating from their strong sexual dimorphism.Phintella is a species-rich genus in Salticidae, and 76 named species with cosmopolitan distribution are included in this genus (World Spider Catalog, 2024).The male-female combination is unknown for about 1/3 number of the species, i.e., 15 species have male-only specimens, and 10 species have female-only specimens.For example, Żabka (1985) described Phintella aequipeiformis Żabka, 1985, based on a male, and P. lucai Żabka, 1985, was based on a female in his monograph of Vietnamese Salticidae.Since then, the opposite sex of each nominal species has remained unknown.Phung et al. (2016) revealed the conspecificity of the two named species based on the results of DNA barcoding and field collection data and synonymized them under the name P. aequipeiformis.
As a case study to shed light on two major problems in our precise understanding of the local species diversity of jumping spiders, i.e., the presence of cryptic species as well as taxonomic confusion caused by strong sexual dimorphism, our present study aims to conduct the species-level delimitation of Phintella and Phintella-like species intensively collected from Vietnam using molecular-based species delimitation using three gene fragments (mitochondrial COI, 16S-ND1, and nuclear 28S) and biogeographical considerations.We used different algorithms to delineate the species boundaries using molecular data.Then, we identified the opposite sex for the putative species for which only one sex has been known.In addition, the divergence time and the diversification rate shift of the biogeographic histories of the major lineages and putative species were also estimated to reveal the governing processes of diversification in Phintella and their phylogenetically related species.

F I G U R E 1
Maximum credibility phylogenetic tree (MCC tree) reconstructed using BEAST and the summary of species delimitation results (part I).The MCC tree reconstructed via BEAST showed that the 22 morphospecies (the PS column in pink) corresponded to monophyletic clades, indicating they are potentially valid species.The species delimitation results of BFD, bGMYC, and ASAP supported the species boundaries of the morphospecies in different degrees.The result of ASAP supported a 30-species scenario, and the result of bGMYC supported a 37-species scenario.For more conservative methods, BFD supported a nine-species scenario, and bGMYC supported a 11-species scenario.
A total of 16 specimens belonging to 14 species of eight genera of Salticidae were used as outgroups (Table S1).Three of these species, Menemerus bivitattus (Dufour, 1831), Siler severus (Simon, 1901), and Telamonia festiva (Thorell, 1887), were collected from Vietnam in this study.All sequences obtained in this study were deposited in GenBank (Table S1).

| Sequence dataset preparation
We align the sequences using the Geneious aligner algorithm in  S1).Such sequences were used for Assemble Species by Automatic Partitioning (ASAP) analysis to test the species delimitation hypothesis based on the two previous datasets (Table S2).

| Phylogenetic analyses and divergence time estimation
Using the concatenated dataset, we used the Bayesian inference to reconstruct the maximum clade credibility (MCC) tree in BEAST 1.10.4(Drummond et al., 2012).We determined each nucleotide substitution model for each gene region (COI, 16S-ND1, and 28S) in jModelTest (Posada, 2008).The input files for BEAST were generated by BEAUti 1.10.4(Drummond et al., 2012).We set the nucleotide substitution model of each gene.We set the substitution rates and clock models of each gene as unlinked, the clock model as an uncorrelated relaxed clock with gamma-relaxed distribution, and a speciation model as the Yule Process Tree Model.We chose Yule Process because we followed the protocol from Bodner and Maddison (2012).Other priors of the parameters were set as default, whereas the molecular clock of the time calibration scheme was set as follows: the MCMC chain length was 1.2 × 10 8 with a sampling frequency of 1 × 10 4 .We visualized the results of independent runs in Tracer 1.7.1 (Rambaut et al., 2018) to diagnose the convergence of independent runs and examine the sufficiency of effective sample sizes (ESS > 200).We performed a burned-in of 40% of trees to generate the final MCC tree using TreeAnnotator 1.10.4(Drummond et al., 2012), which was displayed by using Figtree version 1.4.4 (Rambaut, 2018).
We used these three-time calibrations to estimate the divergence time of the phylogenetic tree.

| Species delimitation analyses using molecular data
We used three approaches for species delimitation, including one barcode gene-based method, ASAP (Puillandre et al., 2020), one Bayesian implementation of the General Mixed Yule Coalescent analysis (bGMYC, Reid & Carstens, 2012), and one Bayesian coalescent-based method, Bayesian Phylogenetics and Phylogeography (BPP, Yang, 2015).We used the concatenated dataset for BPP and bGMYC analyses, and only the COI matrix was used for ASAP.
The bGMYC is the Bayesian version of GMYC (Pons et al., 2006).
The GMYC uses only one ultrametric gene tree as an input and searches the threshold at which branching patterns represent coalescent or speciation events (Emmanuel et al., 2015;Pons et al., 2006).The bGMYC provides the means to use the posterior distribution of multiple ultrametric trees instead of a single tree (Emmanuel et al., 2015;Reid & Carstens, 2012).The 1201 trees obtained from post-burn-in trees generated from BEAST analysis using LogCombiner 1.10.4(Drummond et al., 2012) were used to perform bGMYC analyses.The bGMYC was conducted in the R program using the package "bGMYC."We ran bGMYC under each of the 1000 tree topologies, with 5 × 10 4 MCMC generations, burn-in of 4 × 10 4 , and thinning per 100 generations.ASAP analysis applies a hierarchical clustering algorithm that only uses pairwise genetic distances for building partitions from a DNA sequence dataset of a group of individuals.We performed ASAP at the web interface (https:// bioin fo.mnhn.fr/ abi/ public/ asap/ ) by selecting the JC69 model for the COI dataset.We set the default parameters in "Advanced Option." The partition with the lowest "ASAP score" was selected as the best under a given substitution model.
We used the Bayesian coalescent-based method, BPP (Yang, 2015), to delimitate species boundaries based on the BEAST tree topology.The species status was determined based on the monophyly of the clades under similar morphology and geographic distribution.Accordingly, nine monophyletic groups were set as putative species, and the species of each monophyletic clade was estimated using BPP v 4.4.0 (Flouri et al., 2018).We assumed a guide tree topology from the BEAST tree results and then applied the gamma prior for theta and root tau (gamma 2 1000).BPP analysis was processed using an MCMC chain of 2 × 10 5 after discarding the first 8000 steps.

| Diversification rate shift and ancestral range reconstruction
Based on the species delimitation results, we prepared a final MCC tree from BEAST analysis of the proposed species.We first kept one tip but dropped the other tips in a putative species (species estimated via ASAP and bGMYC as the most specious estimation) and kept the closely related outgroup in the tree.Then, diversification rate shift estimation was conducted using BAMM-2.5.0 (Rabosky, 2014) with four chains, 10 million generations, and simulatePriorShifts = 0.
Priors were set as the recommendation of BAMMtools (Rabosky et al., 2014) for the final MCC tree.Finally, visualized the results by R package BAMMtools (Rabosky et al., 2014) to identify the rate shifts in the MCC tree generated via BEAST using the lineage-through-time (LTT) approach (Nee et al., 1994).LTT estimations were conducted using the R package ape (v5.7.1, Paradis & Schliep, 2019) for the MCC tree and 1000 randomly selected trees from BEAST analysis, then all LTT plots and evolutionary rates visualized by BAMMtools were incorporated by using ggplot2 (Wickham, 2011).The effective sample size (ESS) was evaluated by R package CODA (Plummer et al., 2006); furthermore, functions from BAMMtools were used in the analysis below: we determined the expected number of shifts by computeBayesFactors() and summary.bammdata(), and then the credible set of diversification rate shift configurations were estimated by credibleShiftSet() with the expected number of shifts set to 0. Finally, we visualized the speciation rate by plot.bamm().The ancestral ranges of each node in the ingroup were also reconstructed using the currently known distribution of the putative species (Table S1).We conducted ancestral range model selection to identify the best statistics for reconstruction.The possible models are DIVALIKE, DEC, and BAYAREALIKE implemented in BioGeoBEARs with a j-parameter representing the long-distance dispersal scenario (Matzke, 2013a(Matzke, , 2013b(Matzke, , 2014)).When the best model was identified, ancestral range reconstruction was conducted under BioGeoBEARS methods in Reconstruct Ancestral State in Phylogenies (RASP V.4, Yu et al., 2015).The geographic divisions used in the analysis were based on Bain and Hurley (2011).We define our geographic areas in the RASP analysis according to the number of species included in our analyses and the approximate geographic areas.We divided Vietnam into Northern, Central, and Southern Vietnam.For the areas outside of Vietnam, we designate India, Pakistan, and Sri Lanka as one area; and the Philippines and Indonesia as one area.We set Hawaii, East Asia, and Africa, respectively, as one area.

| Divergence time and diversification rate shift estimations
By using the MCC tree from BEAST analyses and keeping one tip per putative species (37 species scenarios in bGMYC), we show the estimated divergence time of each clade and their inferred diversification rate for a given period, which we did not detect a significant rate shift in our ingroup; rather, the diversification rate is higher during mid-Miocene, ~15-16 MYA, then slowly decline (Figure 3b).Our Other diversification rate shift configurations have very low posterior probabilities (Figure S1), and thus were not considered in our study.This result could be biased because the species are collected from Vietnam.However, the reconstructed ancestral area indicated a SE Asian continental origin (and its neighboring island systems) scenario of the putative species recognized in this study.

| Ancestral range and biogeographic process reconstructions
By contrast, the "non-Phintella" clade (Clade 9) has been primarily diversified in Northern Vietnam and the Philippines + Indonesia (AE).
Moreover, the plot of biogeographic events showed that the governing biogeographic process for the diversification of our ingroup species was dispersal.This result indicated species diversification was due to dispersal and expansion to the adjacent geographic ranges because the long-distance dispersal parameter (J) is insignificant.

| Species delimitation hypothesis and taxonomic remarks
The BFD species delimitation hypothesis with nine species is overly conservative and unreasonable because it proposes the merging of multiple morphological species that can be distinguished based on important characteristics found in somatic and genital morphology.
Apart from the BFD hypothesis, 22 ingroup putative species (PS1-PS22) collected from Vietnam are recognized by considering the agreement between ASAP and bGMYC delimitation (Figures 1 and   2).A total of 11 species are considered to be undescribed species.Of these, the conspecific male-female combination was elucidated for seven species (PS8, PS10, PS11, PS13, PS18, PS19, and PS20), while the only male was found for two species (PS5 and PS16), and the only female was found for two species (PS3 and PS14).Our divergence time estimation showed that the major divergence occurred in the middle to late Miocene but with our significant diversification rate shift.Pleistocene events have little effect on species diversification in our ingroup species (Figure 3).The biogeographic process that increases the number of species through time in our ingroup was the range expansion events via dispersing and diversification to the adjacent ranges in SE Asian landmasses (Figure 4).
The ingroup species were assigned to Clade 3 and Clade 9. Clade 3 consists of 18 putative species from Vietnam, including eight named species and two cases of the merging of named species, namely, PS22 (Phintella bifurcilinea + P. debilis), PS21 (P.monteithi), PS12 (P.aequipeiformis), PS17 (P.vittata + P. suavis), PS15 (P.cavaleriei), PS9 (P.liui), PS7 (P.sancha), and PS6 (P.lepidus); the remaining 10 species in Clade 3 are likely new species (Figures 1 and 2).Considering that P. bifurcilinea is the type species of the genus Phintella, the members of Clade 3 can be reasonably assigned to Phintella.The male and female specimens of Phintella bifurcilinea (Bösenberg & Strand, 1906) and those of P. debilis (Thorell, 1891) were partly mixed to form a single clade (PS22) with a high node support and a relatively long basal branch in ASAP and bGMYC.In addition, the male and female specimens of Phintella vittata and the male specimens of P. suavis were recovered as a single species (PS17).Phintella vittata (Koch, 1846) and P. suavis (Simon, 1885) were described based on different sexes, i.e., the former based on the female holotype from Bintang (Java), Indonesia, and the latter based on the male holotype from Malacca, Malaysia, thereby hindering us to discuss the conspecificity based on the direct comparison of the two species.Therefore, accepting the DNA-based merging of the two named species is reasonable.In Clade 9, multiple species (including PS2, PS3, and PS4) of Clade 11 were morphologically assigned to Phintelloides.
Among the 22 putative species in the ingroup taxa, we identified the opposite sex of two species: Phintella monteithi Żabka, 2012 (male-based species), and Lechia squamata Żabka, 1985 (femalebased species) were reasonably inferred based on molecular phylogenetic analyses (Figure 1).Candidate new species were also recognized from Vietnam (Figures 1 and 2).observed in PS18 and PS19, where both species co-occurred in Pu Mat (Table 2).Although we have newly recorded the morphology of males and females in nine species, we additionally recorded four species in that we only collected one sex.PS3, which is closely related to P. versicolor, and PS14, which is closely related to P. cavaleriei, have only one female specimen in each species.PS5, which is closely related to P. lepidus, and PS16, which is closely related to P. vittata, have only one male specimen in our collection.PS20 was as the sister species of P. arennicolor.
The former phase corresponds to the period after the Middle Miocene Climate Transition, a global drastic change from a warm/wet phase to a cold/dry phase (Methner et al., 2020).This climate change likely modified the spatial distribution of vegetation types (Wang et al., 2019) and their linkages and temporally formed several suitable habitat patches temporally during the transitional period (such as open/sparse forests, forest edges, and bushes).No long-distance dispersal was inferred for the present dataset (i.e., J-parameter = 0), which indicated that dispersal across different landmasses may not always be the case.Rather, regional dispersal and vicariance events could occur at relatively small geographic scales, although jumping spiders could conduct "ballooning" (Jiménez-Valverde et al., 2010).

F
I G U R E 3 Divergence time estimations and diversification rate under the lineage through time approach using the species delimitated via ASAP.(a) The divergence time of the Phintella and the ingroup started at 18.36 MYA.The major clades in Vietnam started diverging after 15.52 MYA.(b)The diversification rate, speciation rate, and extinction rate are demonstrated in the lineage through time plot, the gray lines represent another 1000 LTT plots from BEAST analysis, and shades of each color represent confidence of evolutionary rate calculated by BAMMtools, i.e., the speciation events were not changed through the geological time from Pliocene to Pleistocene (5.00-2.50MYA).(c) The speciation rate was calculated via BAMM.There is no shift in the best diversification rate shift configuration (posterior probability = 0.8).
PS8 and PS9 are sister allopatric species, which are collected from Chu Yang Sin National Park and Vu Quang National Park, respectively.PS10 and PS11 are partly and allopatrically distributed sister species as they cooccurred in the Phu Quoc National Park.The same pattern was F I G U R E 4 The ancestral range reconstructions and the estimation of historical biogeographic events.The RASP analyses show that the best explanatory model is the DEC model (AICc wt = 0.60).The ancestral range reconstruction results under the DEC model show that longdistance dispersal (the J-parameter = 0) did not govern the diversification of the species, instead, short-range dispersal was controlling the species diversification in our target taxa.
were restricted to a narrow geographic range.These facts indicate that Phintella and Phintella-like spiders' species diversity in Vietnam may still be underestimated.The putative species recovered by multiple species delimitation analyses based on the dataset of the three genes (COI, 16S, and 28S) was recovered well by ASAP analysis based on the COI dataset alone (Figures1 and 2).Similarly, COI-based DNA barcoding (or "mini-barcoding" targeting a shorter COI gene marker:Meusnier et al., 2008;Vamos et al., 2017) can also be useful to preliminarily elucidate undescribed species and unknown conspecific male-female pairs in Phintella and its related genera.Apart from the Vietnamese species, additional species other than our targeted geographic areas, for example, China, the Philippines, and other Southeast Asian countries, should be included in future studies to illustrate the evolution and diversification histories of Phintella species.AUTH O R CO NTR I B UTI O N SLuongThi Hong Phung: Conceptualization (equal); data curation (equal); funding acquisition (equal); methodology (equal); resources (equal); writing -original draft (equal); writing -review and editing (equal).Yong-Chao Su: Data curation (equal); formal analysis (equal); funding acquisition (equal); methodology (equal); resources (equal); supervision (equal); visualization (equal); writing -review and editing (equal).Takeshi Yamasaki: Data curation (supporting); investigation (supporting); methodology (supporting); resources (supporting); writing -review and editing (equal).Yi-Yen Li: Formal analysis (equal); visualization (equal).Katsuyuki Eguchi: Conceptualization (equal); funding acquisition (equal); methodology (equal); project administration (equal); resources (equal); supervision (equal); writingreview and editing (equal).